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Abstract. A probability model is presented for the dynamics of mutation- 
selection balance in a haploid infinite-population infinite-sites setting suffi- 
ciently general to cover mutation-driven changes in full age-specific demo- 
graphic schedules. The model accommodates epistatic as well as additive 
selective costs. Closed form characterizations are obtained for solutions in 
finite time, along with proofs of convergence to stationary distributions and a 
proof of the uniqueness of solutions in a restricted case. Examples are given 
of applications to the biodemography of aging. 



1. Introduction 

Arguments from the mathematical genetics of mutation-selection balance figure 
broadly in evolutionary theories of senescence. Available formal models, however, 
do not cover cases brought to the fore by recent progress in biodemography £Q. In 
this paper, we present a rigorous general model encompassing these cases, prove 
results concerning existence, uniqueness, and convergence, obtain closed-form rep- 
resentations for solutions to the model, and give examples of its application to 
questions in the demography of aging. 

The whole mathematical theory of natural selection may be divided into three 
parts: positive mutations, neutral mutations, and deleterious mutations. Positive 
mutations may be thought to add up to an optimal adaptation, at least under 
some conditions, and they are generally studied in that context by demographers. 
Neutral mutations have their primary effects in alleles which drift randomly to 
fixation. Deleterious mutations, the focal subject for theories of aging and for 
this paper, are expected never to achieve fixation in populations, except, through 
founder effects, in very small populations. Their influence in large populations 
derives from their persistent reintroduction and slow meander to extinction. 

Sir Peter Medawar 0, in 1952, descried an explanation for senescence in the 
accumulation of deleterious alleles with age-specific effects, given the declining force 
of natural selection with adult age. W. D. Hamilton |3j presented expressions 
for this declining age-specific force, helping others quantify the resulting balance 
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between mutation and selection. B. Charlesworth 0] analyzed the dynamics of 
age-specific selection. His work guides the thinking of many experimentalists. 

At stake are the cumulative effects of numerous mildly deleterious mutations 
showing up at some large collection of loci. In our setting, the genotypes determine 
full age-specific schedules of mortality and fertility, and the effects of a mutation 
have to be represented as a perturbation of a whole function of age. A rigorous 
treatment demands that mutations correspond to points in abstract spaces, such 
as function spaces. Relationships between our work and the large literature on 
mutation and selection reviewed by Burger are discussed in Section |H1 

Up to now, researchers have relied on linear approximations to cost functions 
and restricted their representations of the age-specific effects of mutations to styl- 
ized patterns like step-functions. Intriguing results have been obtained. Some are 
discussed in Section The linear analysis, however, can be deceptive, and the 
stylized patterns are remote from realistic portrayals of gene action. Cases chosen 
for analytic tractability give a misleading picture of the full range of possibilities. 

Our model is an infinite-population, multiple-sites or infinite-sites model in con- 
tinuous time. The dynamical equation is a fairly standard one, but the space of 
mathematical objects to which it applies is novel. Our model allows a highly flex- 
ible specification of pleiotropic gene action. It is especially suited to demographic 
applications with mutant alleles affecting age-specific schedules. The model is a 
haploid infinite-population model with no recombination. A parallel model with 
free recombination, introduced in Section OH will be developed in a future paper. 

Our contribution is to allow large numbers of interacting genes to make small con- 
tributions to a continuum of linked traits. Traditional analyses which recognize in- 
dividual alleles (thus admitting, in principle, arbitrary configurations of pleiotropy) 
are amenable only to small numbers of loci; quantitative genetics, which reduces 
the contributions of individual genes to a continuum, reduces the complexity of 
pleiotropy to covariance matrices. 

Although multi-locus models without recombination like our own can be formally 
imbedded in single-locus models, this imbedding will not generally yield useful re- 
sults. When a multilocus model is translated into the single-locus framework, it 
brings along an extra structure of transition rates, whose complexity grows ex- 
ponentially with the number of loci. When the number of loci is large or, as in 
our model, effectively infinite, this extra structure overwhelms the single-locus in- 
frastructure. In our function-space setting, the formal embedding itself also poses 
difficulties. As a consequence, results for single-locus models are mainly helpful as 
analogies. 
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Unlike most models of which we are aware, our model comfortably accommodates 
epistasis. (A very different approach to epistasis, in the two-allele setting, may be 
found in UJ.) The selective cost of a mutant allele can depend on the configuration 
of other mutant alleles present in a genome. This property is critical to the study 
of senescence, even without special assumptions about interactions among genes, 
because the fitness costs of cumulative demographic changes are not linear. 

We are able to obtain closed-form representations of the entire time path of 
solutions to our dynamical equation fTheorem l3.1|l . Our results are not restricted, 
like much previous work, to limiting states and equilibrium distributions. We give 
proofs of convergence over time fTheorem I4.1JI . and set machinery into place to 
compute rates of convergence and to cope with changing fitness conditions as well. 
In Section |5j we present some results about the asymptotic behavior of solutions. 
Theorem 15.11 gives sufficient conditions for the numbers of certain classes of mutant 
alleles to increase without limit, generalizing the well-known "error threshold" (cf., 
0.) In Sectional we derive the Poisson limit for the non-epistatic case, as well as 
proving uniqueness of the solution. In the general epistatic case we do not yet have 
a proof of uniqueness. In Section we discuss some implications of our results for 
the theory of longevity. In Section |S1 we review earlier work on related problems. 



2. The model 

We consider an infinite population subject to mutation and selection. There is 
a complete, separable metric space M of potential mutations, on which is defined a 
boundedly finite Borel measure v. (In other words, v assigns finite mass to bounded 
sets; together with the assumptions on M, this condition implies that v is cr-finite.) 
We refer to this measure as the "mutation rate"; for any set B, the quantity v(B) 
represents the rate at which there spontaneously arises a mutant allele from B. Our 
picture is one in which new mutant alleles are steadily arising, each one tagged by a 
corresponding point of M. For convenience, we identify the tag with a description 
of the effects that the mutant allele produces: for instance, a function on the non- 
negative real line R + giving the increases in mortality attributed to the action of 
that allele at each age. 

The space of "genotypes" S is identified with the integer-valued boundedly finite 
Borel measures on M, with a topology to be described shortly. An element of S has 
the form S mi , where the rrii £ M are not necessarily distinct and the number of 
rrii in any bounded subset of M is finite. The notation 8 X stands in general for a unit 
mass at the point x in the space to which x belongs. Each genotype represents a set 
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of mutant alleles that an individual may carry. The "null genotype" has wild-type 
alleles at every locus and carries none of these mutant alleles. 

The state of the population at time t is denoted Pt , which is a Borel probability 
measure on the measures in 9- Thus Pt is the distribution of a random measure 
[El El- The evolution of the population is presumed to be so slow that it can be rep- 
resented as occurring in continuous time, without reference to discrete generations. 

To each genotype g we assign a "selection cost" S(g); S is a continuous function 
from 9 to R + . (Including negative costs would be feasible for the finite-time so- 
lutions, at the expense of slightly more complicated statements for theorems.) In 
applications, costs will typically be decrements to growth rates, in effect measuring 
fitness on a logarithmic scale. 

We normalize costs so that S vanishes on the null genotype, and vanishes for 
no other g. On M we write S(m) for the cost of the singleton g = 5 m . When S 
is linear, so that S(g + S m ) — S(g) is independent of g, the model is additive, or 
nonepistatic. 

Any measure P on 9, like P t , may be determined by the expectation values it 
assigns to a suitably rich collection of functions F from 9 to K, as specified below. 
For brevity, we write PF or P(F) for the expectation value PF — J g F(g)dP(g) 
of any measurable function from 5 to 1 such that J \F(g)\dP(g) < oo. Since 
genotypes are measures, we can also write gf = g(f) = J f (m)dg{m) = ^ f(rrii) 
when / : M — > M, and g = ^ 6 mi . 

Our dynamic equation for P t is 



The meaning of the equation is readily described when F is the indicator function 
of a set G of genotypes. The first term inside the integral measures the rate at which 
the population is flowing into the states in G out of all sorts of other states because 
of the addition of a mutant m that lines up just right to enter G. The second term 
inside the integral measures the rate at which population flows out of G because 
of new mutations. The remaining two terms measure the effect of selection. The 
proportional rate of change in mass of the population in G equals the difference 
between the average fitness cost of genotypes in G and the average fitness cost of 
the whole population. 

Measuring fitness relative to the changing average fitness of the whole popula- 
tion keeps total mass constant and lets the measure Pt represent the probability of 
finding a randomly selected individual in a given state, modeling population distri- 
bution rather than population size. While our equation may be novel, it is strongly 
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analogous to standard mutation-selection dynamics on quantitative traits, such as 
those given in equation V. 2. 11 of 5 . 

When mutations are identical (so that M comprises only a single point) we have 
the "mutation-counting model" going back to Kimura and Maruyama ^U|j whose 
history will be described in Section|Hl A genotype is specified by a natural number, 
the number of mutant alleles present in it, and Q becomes 



In the non-epistatic case, where S is additive, the mutation-counting model (or its 
discrete-time counterpart) has a Poisson distribution with parameter i//S(l) as its 
stationary distribution. 

For general S, the counterpart of the Poisson distribution is a Poisson ran- 
dom measure. For the non-epistatic case, Theorem 16.11 establishes conditions for 
uniqueness and convergence to a stationary distribution given by a Poisson random 
measure with intensity (1/ ' S(m))dv(m), the measure on M whose Radon-Nikodym 
derivative with respect to v is 1/S. (A Poisson random measure assigns a Poisson- 
distributed random integer mass to each measurable set, the mean of the mass 
assigned to a set is the intensity measure of the set, and the random masses of 
disjoint sets are independent random variables.) The general theory 0[ Chapter 7] 
takes care of technical details. Even in the non-epistatic case, only rather special 
starting states lead to the Poisson limit. In the epistatic case, covered by Theorem 
13.11 asymptotic distributions, when they exist, may not be Poisson. 

We need a suitable notion of weak convergence for boundedly finite random 
measures. Following Appendix A. 2 of UJ, we equip the space S with the metriz- 
able w-topology under which a sequence of measures <?i, (72, ■ • • € S converges to 
a measure g G S if and only if lim„ g n (f) = s(/) f° r each bounded continuous 
function / : M — > R that is supported on a bounded set. A sequence of probability 
measures Qi,Q2, ... on S (that is, a sequence of distributions of boundedly finite 
random measures on M) converges weakly to the probability measure Q on S with 
respect to the ^-topology if and only if \\m n Q n (F) = Q(F) for every bounded 
w-continuous funtion f on S. This turns out to be equivalent to the requirement 
that lim„ Q n (F) = Q(F) for all F of the form F(g) = e~ 9 ^ for some continuous 
boundedly supported / : M — * K + . This class ? is the sufficiently rich class of 
functions required for our expectation values to determine our measures: Equality 
of expectations for F in J implies equality of expectations for all bounded Borel- 
measurable F [HI Section 6.4]. 



(2) 



dP t {n) 
dt 



vP t (n - 1) - vPtin) - P t (n) S(n) - ^ S{m)P t {m) 



m 



) 
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We must bear in mind that ordinary theorems guaranteeing existence and unique- 
ness of solutions to differential equations do not extend to our abstract setting. The 
derivatives on the left-hand side of Q might not exist, and, in the presence of un- 
bounded S 1 , we might have infinity minus infinity on the right-hand side. Our proofs 
are constructive, so the meaningfulness of the equation will follow from the prop- 
erties of the proferred solutions. We prove a reasonable version of uniqueness in 
the non-epistatic case. In the general epistatic case, we have not yet ruled out mul- 
tiple alternative meaningful solutions; there could be a complicated mathematical 
question lurking here. 

3. Existence of solutions 

We express the solution in terms of a certain random meaure on M. x R + . We 
let IT denote the Poisson random measure on this space with intensity measure 
v ® Lebesgue. Define a time-homogeneous S- valued Markov process (X t )t>o by 

(3) X t :=X + [ S m dU(m,u), 

;j»tx[o,i] 

where Xq is a random measure with distribution Pq, independent of II. Each real- 
ization of X t may be pictured as a discrete set of points, possibly with duplication; 
as time passes, new points accrete. The cost function S could be allowed to depend 
on time, but we keep time-independent notation here. 

Theorem 3.1. Suppose that there is a positive T such that 

(4) Ecxp (- J S(X u )dv\ S(X t ) < oo 

for all t € [0, T). Then the equations Q have a solution on [0, T), given by 



(5) P t F 



Proof. Define a linear operator on the continuous functions on the genotype space 

by 

(6) AF = J [F(- + S m ) - F(-)] dv{m) - S(-)F(-). 

Given an integrable function u(t), put P t = cxp (^— J Q * a(u) duj Pt- If we can ar- 
range for a(t) to equal the average selective cost PtS, then, thanks to the chain 
rule, the derivative of P t F must equal Pf AF. 

The operator A may be unbounded if v has infinite total mass, but it is well- 
defined on the class J and is the generator of a sub-Markovian semigroup {T t )t>o- 



E 


exp J* * S(X U ) du^ 


F(X t ) 


E 


exp(- J*S(X U ) 


du^j 





MUTATION-SELECTION MODEL 



7 



By the Feynman-Kac formula Section III. 19], T t may be described as 



(7) T t F(g) = E 



cxp (- J S(g + X U - X ) du^j F(g + X t - X ) 



Now, the semigroup (T t )t>o solves the forward equation ^T f F = T t (AF) and 
it follows that PtF — PoT t F, which equals the numerator of (J5J. By the condition 
on T, P t S is finite on [0, T), equalling the derivative of PtI, so that we may put 
P t F = P t F/P t I and achieve <r(t) = P t S < oo on [0, T). □ 

4. Representations 

Although our solution J5J may look abstract, as long as f(M) is finite PtF can 
be expressed as a series expansion whose terms can be evaluated by multiple inte- 
gration. We now derive this expansion, which makes direct calculations feasible in 
applications. When i/(M) is finite, we order the points put down by II according 
to their arrival times t(1),t(2) . . . and write Y„ := X T i n y Let J n be the indica- 
tor function of genotypes with exactly n (possibly overlapping) points: J n (g) = 1 
if g(M.) — n, and otherwise. Our series expansion for P t F will take the form 
^2 n PtJ n F. We write x Ay for the lesser of any two quantities x and y. Re- 
newal theory calculations turn j£| into a handy formula for probabilities of rt-point 
genotypes: 

Theorem 4.1. Suppose v{JA) < oo and Pq puts unit mass at the null state 0, with 
S(0) = 0. Then the solution J5J may be written as P t F = P t F j P t \, with 

(8) P t J n F = u(M) n e- v( - M ^E [(S(Y 1 ) . . . S^y 1 H t , n F(Y n ) 

Here H t n is a conditional probability defined in terms of independent unit-rate 
exponential variables Z 2 , ■ ■ ■ by the formula 

(9) H t<n - P Z i/ S ( Y i) <t\Y ll ...Y n } 

If ^2 ^(M)"E[((S'(Yi) . . . S(Y n ))^ 1 ] is finite, P t converges in distribution as t goes 
to infinity. If the sum is infinite, PtJ n goes to zero for all n. 

Proof. Consider the numerator of JHJ with J n F in place of F. The integral inside 
the exponential is the sum of terms S(Yj)(r(j + 1) At — r(J)) for j from l to n. The 
factor J n restricts the domain to the event {X t (M) = n}, an event with probability 
e~ vt v n t n /nl, where we write v for i/(M). Conditional on this event, the Y's are 
independent of the r's, and the r's are distributed like the order statistics of a 
sample of n uniform random variables u\ . . . u n on [0, t] which may occur in any of 
nl orderings. Put it n +i = t. 
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To obtain the expectation over the r's, we evaluate the integral 

(10) nlr n / ••■/ exp I - S(Yj)(u i+ i - m) > du„ ■ ■ ■ dmdui, 

JO All J«„_i { i=1 J 

The change of variables zi = — ui transforms this integral into the product 
n\/{t n S{Y 1 )...S(Y n )) times 

(11) J - J (siY^e- 3 ^) . . . (S(y n )e- S ^") dz x dz 2 . . . dz n , 

The integrations range over all non-negative z\ . . . z n such that z\ + . . . + z n < t, 
yielding the exponential probability expression H t , n - Closed- form formulas for H 
are given in [TJ Ch. 1, 13.12]. The probability e-"v n t n /n\ times n\/{t n S{Y{) . ..S{Y n )) 
times H tt n gives ©. 

We bound PtJ n S by isPtJ n -i, noting that Ht n < Ht, n -i- Summing over n, 
we find P t S < vPtl < v, verifying the supremum condition for all finite T. The 
factors of e~ ut in the numerator and denominator of PtJ n F cancel. The conditional 
probability H t , n , is monotone increasing in t toward a limit of 1 for each choice of 
n and Y\ , . . . , Y n . Hence the limit claim follows by monotone convergence. □ 

In demographic applications we are typically interested in counting the average 
number of mutant alleles of a given type that a randomly chosen individual would 
bear. For B a measurable subset of M, write Rt(B>) for the expected number of 
mutations from B at time t; that is, Rt(B) = J g g(B)dPt(g). For special starting 
states, we can obtain a closed-form density for R t . 

Theorem 4.2. Suppose the starting distribution P$ is a Poisson measure with 
intensity ttq. Then the measure Rt has the form Q(m)dv{rn) + T]t(m)d7ro(m) where 



E[cxp(- /p S(X u )du) /„' exp(- j'[S(X u + S m ) - S(X u )]du)dT] 
1 m ~ E[exp(- ^ S(X u )du)} 

, , E[exp(- f*S(X u + 6 m )du)] 

Tit ( 777 ) = — — — 

E[exp(- $lS{X u )du)] 

Proof. When the initial distribution is Poisson, the entire process X t , including 
Xq, is defined from a Poisson random measure £ = II + (X ,Sq) on the product 
space M x R + with intensity measure H = v ® Lebesgue + 7r €5 <5o- The local 
Palm distribution for the Poisson random measure £ at (to, t) in JVC x M + is the 
distribution of £ itself augmented by an atom at (to, t) Example 12.1(b)]. For 
any non-negative bounded Borel-measurable function G(m,r, £) the Palm integral 
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formula Proposition 12.1. IV] makes 

(13) E [ G{m,T,£)d£(m,T) = [ EG(m,T,( + 8 {mjT) )dH{m,T) 



The integrals are taken over M x R + , and E operates on £. Fix t and B C M and 
choose the function G to be 

(14) G(m, t, £) = exp f- ^ £(XJ du^ l {meB} l {r < t} . 

Bear in mind that X u is a function of £. With this G, plugging into Equation JSJ, 
PtXt(B) is given by the left-hand side of 113(1 . On the right-hand side, the extra 
atom at (to, t) changes the argument of the exponential function inside l(14|) into 
— J Q T S(X u )du — J T S(X U + 5 m )du. The first term in H , which is z/ ® Lebesgue, calls 
for integration over r and gives the contribution in the numerator of £ in (|12fl with 
respect to z/. The second term in H puts r equal to zero and gives the contribution 
in the numerator of r/ with respect to 7To. The denominator in £ and ry is a constant 
independent of the set B. It converts Pt to Pt- The indicator function in G arranges 
that the measure Rt{B) is obtained by integrating over B, so £ and rj are indeed 
Radon-Nikodym derivatives for R t as claimed. □ 

When the process Pt starts from the null genotype we set tto = 0. Equations 
(|T2"f allow us to compare the influences of different cost functions: 

Corollary 4.3. Assume that the conditions of Theorems \3.1\ and \4-S\ are satisfied, 
and suppose that S is sub-additive; that is, S(g + g') < S(g) + S(g'). Define 
the corresponding additive cost function S(g) := / S(S m )dg(m). Let Pt and Pt be 
corresponding genotype distributions produced by JSJ. Then P t F > P t F for any 
linear F of the form F(g) = g{f), where f is nonnegative, measurable, and has 
bounded support. 

Proof. The sublinearity of S and the linearity of S imply 

E[exp(- /* S(X u )du) Jexp(-(t - T )S{6 m ))dr] 
Ct(m) > — — -. — - 

(15) ' Ecxp(- fiS(X u )du) 

I _ e -S{m)t 

= q7 "i = Ct(m). 

b (to) 

Similarly T)t(m) > e~ s ^ m ' t: qt{m). The result follows from the special case of the 
Palm integral formula known as Campbell's Theorem |9., (6.4.11)]. □ 
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5. Asymptotic Behavior 



In contrast to the additive case, genotypes subject to subadditive cost functions 
may tend to explode. In age-structured models, the total effect of mutant alleles 
acting after some given age is limited, regardless of how many of them may ac- 
cumulate. If the rate at which some class of mutant alleles is generated exceeds 
any countervailing selection-cost increment which they may incur, the number of 
mutant alleles in that class may be expected to grow without limit. 

Theorem 5.1. Assume the conditions of Theorem are satisfied. Let B C M 
be a subset with finite v-mass. Suppose < S(g + g*) — S(g) < s for all g and for 
all those g* with masses only at points in B, that is, with g*(B) = g*(M). Let J* 
be the indicator function of the set of genotypes with g(B) = n. Then s < v(B) 
implies that PtJ„ goes to zero for every n as t goes to infinity. 

Proof. We write our Poisson process X t as X^ + X r t , where X t * is the restriction 
of X t to B and X"[ is the remainder. These components are independent of each 
other. Let U :— inf{u : X t (B) > 0} be the arrival time of the first point in B, an 
exponential random variable with mean \jv[B). We have 



To bound Pt J„ , we write the numerator of j£| , Pt J% , as the expectation of a prod- 
uct of three factors, exp(- J Q S(Xl l )du), exp(- f Q S(X U ) - S(X^)du) and J*(X t *). 
The second factor is bounded above by 1 and the third factor is independent of 
the first. The denominator of JSJ), PtI, has the same first factor, the same second 
factor, bounded below by exp(— s(t — U) AO), and a third factor identically equal 
to 1. Using independence, we may cancel the expectations of the first factors in 
numerator and denominator, so that -P*J* is less than or equal to the quotient 
of EJ*(X t *) and Eexp(— s(t — U) AO). Writing v* for v(B), this quotient equals 
{v*t) n /n\ divided by (^*e^*~ s ^ - s)/{v* - s). The quotient goes to as t — > oo 
for every n. □ 



In the non-epistatic case, when the cost function 5* is additive, a proof of unique- 
ness and an eminently computable formula can be obtained which lead to conditions 
for convergence as t goes to infinity: 

Theorem 6.1. Suppose that S is an additive (nonepistatic) cost function such that 
the expectation value v{S A 1) is finite and suppose that Pq is an initial probability 
measure such that PqS is finite. Then the equations Q have a unique solution on 



(16) 




6. Non-epistatic cost functions 
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[0,oo). A random measure chosen according to Pt may be represented as the sum 
of two independent random measures. The first component is a Poisson random 
measure with intensity (1/ S(m))(l — e~ s ^ m > t )du[rn) . The second is the initial mea- 
sure Pq, tilted by the weighting e~ tgS . That is, the second component Q t satisfies 
Q t F = QtF/Qd with 

(17) Q t F = J e' s ^ t F(g)dP (g) 

If v is finite, this solution is identical with that given in Theorem \S.1\ 

Proof. Linearity of S allows us to transform Equation Q into a first-order linear 
partial differential equation. Suppose we are given an mtegrable non-negative func- 
tion o~(t) which serves as a candidate for PtS. Let z be a positive real number and 
let / on 9 be a bounded nonnegative function with bounded support. We take our 
test functions F now to be of the combined form F(g) = e~ 9 ^~ zS ^- 9 \ We write 
h(t, z) for the real function which will turn out to be log(P t F) satisfying given 
boundary conditions 

(18) h(0,z) =v(z) = log P F. 

Thanks to the form of F and the linearity of S, the expression J (F(g + S m ) — 
F(g))dv(m) from Q now equals F(g)£(z), where 

(19) ((z) := y(F{5 m ) - 1). 

Since / is non-negative, |C( 2 )| is bounded by v{f) + (1 + z)v(S A 1). The first term 
is finite because / is bounded with bounded support. The second term is finite by 
assumption. If exp(h) is to satisfy we need h to satisfy the following partial 
differential equation of the McKendrick type familiar to demographers: 



(20) ^-^ = t(z)+o-(t). 

at oz 

We have shown that the right-hand side is well-defined for all non-negative z and 
t. 

We solve (|20|l uniquely for h by exploiting the method of characteristic curves 
to tranform it into a system of ordinary differential equations. The characteristic 
curve passing through the point (t, z) is the line r 1— » (t, i + z — r) ^3 Section 3.2]. 
Defining h(r) := h(r,t + z-r), we get h!(r) = a(r) + ((t + z - r) for < r < t + z. 
Integrating this equation from to t gives 

pt rt+z 

(21) h(t,z) =r)(t + z)+ a(r)dT+ ((r)dr 

J0 Jz 
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The final term in £ is equal to 



(22) 



-t + 



-f(m)—za(m) 



-/(m)-(*+t)«(m) 



/ 's(m) 



We now set P t P = exp(h(t, z)). Additivity of 5 makes the derivative of PtF 
equal the right-hand side of QJ plus (<r(t) — P t S)(PtF). Also, —PtS is the partial 
derivative of /i with respect to z at z = and / = 0, which is the sum of ^(1 — e ) 
and PoSe~ ts / Poe~ ts . Setting —a(t) equal to this sum is therefore the unique 
choice which makes Pt satisfy QJ. Writing out h and setting z = 0, we recognize 
the Laplace functional of the convolution of probability measures specified in the 
theorem. □ 



The first piece of P t clearly converges to a Poisson random measure as long as v/S 
is boundedly finite. But that is not the complete story of asymptotic behavior. In 
general, the influence of Pq in Q t may persist. In the limit, however, we may apply 
Varadhan's Lemma Theorem III. 13] to show that Q t becomes concentrated on 
the set of genotypes of minimum selective cost. 

Corollary 6.2. Suppose Pq and S satisfy the conditions of Theorem W.1V that the 
support supp Pq is compact, and that S is continuous. Let a = inf{5(<7) : g 6 
SuppPo}. If is any open neighborhood of {g 6 suppp) : S(g) — cr}, then 
lim^oo Qt(0) = 1. In particular, if Po{0} > 0, the tilted measure Qt converges to 
So- 



7. Applications to the theory of longevity 

We outline a few of the many applications to the biodemography of longevity. 
We take the space of potential mutations M to be C[0,oo), the continuous real- 
valued functions on R + , supplied with any of the usual metrics corresponding to 
uniform convergence on bounded intervals |15l Section 1.44]. A mutation measure 
v on this space is the distribution of a stochastic process. We base our selective 
cost function S(g) on Equation 4.9 of Charlesworth 4 , p. 140], taking into account 
more recent discussion [161 p. 930]. The cost function is defined in terms of the 
age-specific survival function l x (g) and the age-specific fertility rate f x {g) specific 
to each genotype, along with a conversion factor T, representing a baseline value 
for the length of a generation, and a rate ro representing a population-wide baseline 
intrinsic rate of natural increase usually set to zero in applications. 

(23) S(g) = (1 - J e-^l x {g)f x (g)dx){l/T) 
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7.1. Gompertz hazards. Charlesworth ^2] has suggested a possible origin for 
Gompertzian (exponentially increasing) hazard rates through a process ol mutation- 
selection balance which fits into our generalized model. Members of a species are 
taken to be subject to a common high background age-independent hazard rate A 
plus age-dependent contributions from mutations. Each mutant allele m may be 
represented as a continuous function m{x) of age added onto the hazard function 
for an individual. Charlesworth's elementary models assume constant fertility at 
all ages above an age b of sexual maturity, forgoing any a priori upper age cutoff. 

The selective cost S(g) for a genotype g from Ij2^(l takes the following form when 
time is measured in generations rather than years: 

(24) S(g) = l-J Aexp(-Aa; + Xb - J J m(a)da dg(m)jdx 

This cost function S is a non-additive epistatic cost function. Following established 
practice, Charlesworth substitutes the additive cost function 

(25) S(g) = J (J ( e -^- b ) A l)m{x)dx^ dg(m). 

This function S is an additive approximation to S. 

We first show that under the same premises as ^7] our model confirms the 
same conclusions. With additive costs as in l)25[l. Theorem 16.11 and Corollary 16.21 
give us sufficient conditions for the distribution of genotypes to converge to the 
Poisson random measure with intensity v / S. It suffices that the starting state put 
positive weight on the null state and have compact support and that v and v j S be 
boundedly finite. Then the average of the hazard rates over genotypes will converge 
to A + f M (m(x)/S(m))dv(m), equivalent to ^3 Eq.4a]. 

It is worth mentioning that this expression for the average of the hazard rates 
is not the equilibrium aggregate hazard rate for the whole population, because the 
heterogeneity mediated by the Poisson distribution implies attrition of higher-risk 
genotypes with advancing age. The Poisson expression for the additive genetic 
variance and covariance also require modification for age-specific attrition. 

Charlesworth focuses on translation families of mutations, which we may write as 
m y {x) = m (x~y) with effects only after an age of onset y. With dv(m y ) — v^dy on 
some [b' , oo), he displays choices for mo which make the average of the hazard rates 
into an exact Gompertz-Makeham function A + i/q cxp(— X(x ~b)) on the support of 
v, and others which approximate Gompertz-Makeham shapes for large x. (These 
shapes do not include heterogeneity corrections.) 

We now observe that our generalized model predicts different qualitative behav- 
ior when the additive approximation of (|24|l by (|25|l is not guaranteed to hold. The 
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additive theory predicts that a "wall of death" with an infinite equilibrium mean 
hazard rate appears at, but not before, the age at which reproduction comes to 
an end |171 p. 60]. Theorem 15.11 implies a more dramatic breakdown. The mean 
hazard function can actually reach infinity at ages at which fertility is still strictly 
positive, if the full epistatic cost function S in l|24(l is kept in place of the additive 
approximation (|25(l . The same is also true, if the bounded cost function S is re- 
placed by an unbounded cost function defined, as in Equation 4.12 of @] p. 141], 
to equal the decrement to the intrinsic rate of natural increase resulting from the 
mutations contributing to each genotype. Contrary to additive theory, the "wall 
of death" is not tied to the end of reproduction but involves a fine-tuned balance 
between mutation rates and tapering costs. 

7.2. Gaussian process mutations. We now apply our model to move beyond 
stylized cases and investigate a wider range of possible specifications for the age- 
specific effects of mutations. The cases considered in 17.11 in which a constant 
background mortality imprints a Gompertzian pattern onto increments to the haz- 
ard function, share the property that every mutant is deleterious at every affected 
age. Is this property essential to the imprinting, or can the age-specific force of se- 
lection readily produce the same kind of outcomes with mutants that mix positive 
and negative effects? 

Our framework allows quite general pleiotropic specifications. A natural start- 
ing point is the case of Gaussian processes. The fitness cost for this brief discus- 
sion will be the additive approximation (|25J) . Suppose that the mutation process 
generates mutations proportionately to a positive-real-parameter Gaussian process 
with expectation a(x) and covariance function c(x,x'), conditioned on fitness cost 
bounded away from 0. That is, if we look at the pattern of age effects in a ran- 
domly chosen mutation, it looks like a realization of this Gaussian process, subject 
to S(m) > s > for some s. The overall rate of mutation is a constant v>q. For 
rigorous treatment, we also need to condition on events which keep the resulting 
hazard functions non-negative and insure the validity of the additive approxima- 
tion, but here we shall assume that the choices of parameters keep misbehavior rare 
enough that it can be neglected. 

The average over genotypes of the hazard function is given by 



numerator and denominator of l|26|l . linear functional of m, are both Gaussian 
random variables, we can describe their joint distribution simply by computing 



(26) 




The denominator S(m) is obtained from (|25|l by integration by parts. Since the 
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their covariance. The conditional expectation of m(x) / S(m) conditional on X(m), 
is obtained via linear regression. 

When a(x) = 0, the conditional expectation turns out to be independent of 
S(m), making the integral in (|26|l independent of the bound s (except for a small 
change in the proportionality constant) and equal to i^Km(jx)S(m)/ Var(£(?n)). 
For processes with zero mean, then, we may treat s as 0. As one example, take 
the mutations to be realizations of Brownian motion, so that c(x, x') = x A x' and 
a(x) = 0. The increment to the mean hazard rate then has the form c\ — C2e~ Xx 
for x > b. 

Can any Gaussian mutation process with zero mean generate approximations to 
Gompertz-Makeham hazard functions? The covariance kernel must satisfy \c(x, y)\ < 
c(x, x)/2 + c(y, y)/2. With a{x) = 0, the incremental mortality is bounded above by 
Ci + C2c(x,x) for constants ci and c-i- The mortality thus cannot be exponentially 
increasing over a long range of ages, unless this exponential increase is built into 
the mutation process itself. 

8. Historical Background 

We discuss in this section the relationship of our model to the existing corpus 
of work on related topics. It was J.B.S. Haldane ^B] who articulated the concept 
of mutation-selection balance as early as 1937. Crow and Kimura ^Hj, Ewens |2()j . 
and Kingman [21] give the foundations of the subject. Burger [23 an d covers the 
present state of the art. These authors put only limited emphasis on age structure; 
Charlesworth 0] and j2H] propounds the age-specific side. 

Infinite population models in which fitness is a function of the number but not 
the identity of mutant genes go back to Kimura and Maruyama ^UJ- They state 
discrete-time and continuous-time dynamic equations for special cases which readily 
suggest the general "mutation-counting model" J2J ■ They obtain some closed-form 
equilibrium distributions. Conditions for convergence to stationary states follow 
from a theorem of Kingman [21] , generalizing theorems of Moran [25] [2U • Burger 
pp. 298-308] traces the subsequent history. Markov-chain versions with stepwise 
mutations of identical deleterious effect have Poisson stationary distributions (see 
e.g. Haigh [22] or Durrett [2H1 P- 137]). The Poisson limit is implicit in estimations 
of equilibrium genetic variance |17j . 

Mutation-counting models in the tradition of Kimura and Maruyama are more 
tractable than the general case considered here because they are a kind of multi- 
locus model that can be subsumed under the theory of single-locus models. The 
count of mutant alleles at different sites can be likened to the integer label on 
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a countable set of alleles at a single site, subject to constraints on the non-zero 
interallelic mutation rates. Models defined by various alternative sets of constraints 
have been studied in some detail. 

The most famous of these single-locus models is Sir John Kingman's "House of 
Cards" (HC) described in [21] and [21] ■ Kingman's infinite-population discrete- 
time model posits a single gene with potentially infinitely many alleles. Alleles 
mutate to new alleles at a constant rate; each new allele has a random fitness, 
given by a probability distribution on [0, 1]. The state of the system is given by 
a distribution of fitnesses on [0, 1], and the dynamics are governed by a standard 
evolution equation. Kingman |29j gives the original proof that the distribution of 
fitnesses for the HC process converges to a limiting distribution. This model has 
many descendants, including the "HC-approximation" |30| for stabilizing selection 
around an intermediate optimum. 

Our model differs from HC and its counterparts in four main ways. Mutant alleles 
in HC have no properties other than fitness. Mutant alleles in our model are tagged 
by an effect represented by a point in a general metric space whose specification 
determines the fitness through the impact on demographic rates. In HC there is 
only a single locus. In our model we are concerned with the heterogeneity of whole 
sets of mutant alleles across a large number of loci within the population. Because 
HC includes only a single locus, it offers no possibility for interactions between the 
fitnesses of different alleles. Our model is open to general epistasis. Finally, HC 
is well-suited to the use of Markovian methods and sample-path analysis, whereas 
our proofs require non-Markovian machinery. 

A highly versatile general formulation of single-locus models has been developed 
by Reinhard Burger chapter IV. 2]. His "general mutation-selection model" 
warrants close comparison with our own. Like us, Burger draws mutants from a 
general space, according to an arbitrary distribution. Burger requires his space to 
be locally compact, whereas we allow any complete, separable metric space so as 
to include mutants identified with continuous functions on [0, oo). 

The substantive difference between Burger's model and our own is in their con- 
trasting views of the genome. Burger focuses on a single locus, with (perhaps) 
infinitely many potential alleles. Each individual's genotype is characterized by a 
single quantity and the population is characterized by a distribution on the muta- 
tion space. We take a more synoptic view, watching the (perhaps) infinitely many 
alleles pop up at (perhaps) infinitely many loci, thereby opening up the represen- 
tation of population heterogeneity. In our model, the only mutation process is the 
conversion of an undifferentiated wild type to a random mutant allele, so we need 
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not introduce transition rates between alleles. However, our flexible treatment of 
heterogeneity means that even the description of the state of the system has to be 
more abstract than is customary in population genetics. 

In a different setting, Del Moral and Miclo [3JJ present results which parallel 
our Theorem 13.11 Their conditions are more general in some respects and more 
restrictive in others. While their concerns are remote from biology, they use the 
terminology of "mutation generators" and "adaptation" in their descriptions. They 
prove that the differential equation model which we analyze can be derived as an 
infinite-population limit of finite nondeterministic Moran models for interacting 
particles. A new book |H2| expands this line of investigation. 

We accompany Del Moral and Miclo, on a road that diverges from the Markov 
modeling, branching processes, branching diffusions, and superprocesses which are 
so important to stochastic population theory [jJ3J. Pioneering work in these areas 
by [3j] 023 On] has been followed by extensive results on particle processes and 
measure- valued diffusions with selection, including |571 15511591 13011211 ELI EC51 EH] . 
Such Markov processes, even when they involve selection, are essentially linear. 
Lineages rise or fall at their own rates, according to their fitnesses, independent 
of the outside population. By contrast, the mutation-selection paradigm on which 
we focus has to be nonlinear, since every lineage has negative fitness. The models 
are saved from trivial degeneracy by a renormalization, conditioning the process 
on long-term survival. This ingredient introduces a quadratic nonlinearity into the 
evolution equation, inasmuch as the entire population contributes to the selective 
pressure on each individual, bringing non-Markovian arguments to the fore. 

9. Prospects 

The generalized model for mutation-selection balance presented here can be ap- 
plied widely to settings where age structure matters. Because the model allows mu- 
tations with a mixture of positive and negative effects, it gives scope to some blend- 
ing of ideas about mutation accumulation with ideas about antagonistic pleiotropy. 
It offers a handle on responses to changing fitness conditions through the finite-time 
solutions, along with machinery for treating epistatic cost functions. 

The Palm formula in Section Ej] facilitates the construction of an alternative 
version of our model which, in contrast to Q , allows for free recombination (FR) . 
In line with ECU], EE)]> and ERJ, we postulate conditions on the relative rates of 
recombination, selection, and mutation which lead, in the continuous-time limit, to 
a process in which P t is always a Poisson random measure on S with some intensity 
measure pt . Our results in Section E^ derived in the absence of recombination, 
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give us the form of an equation for pt in this generalized free-recombination model. 
Differentiating l|12(l at t — leads to a representation for pt of the form 



Here Dp is a measure whose density with respect to p at m is E[S(X p +5 m ) — S(X P )], 
and X p is the Poisson random measure with intensity p. Rigorous development of 
this alternative is reserved for a future paper. 
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